Elastocapillary self-folding: buckling, wrinkling, and collapse of floating filaments 
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When a flexible filament is confined to a fiuid interface, the balance between capillary attraction, 
bending resistance, and tension from an external source can lead to a self-buckling instability. We 
perform an analysis of this instability and provide analytical formulae that compare favorably with 
the results of detailed numerical computations. The stability and long-time dynamics of the filament 
are governed by a single dimensionless elastocapillary number quantifying the ratio between capillary 
to bending stresses. Complex, folded filament configurations such as loops, needles, and racquet 
shapes may be reached at longer times, and long filaments can undergo a cascade of self-folding 
events. 



I. INTRODUCTION 

A common strategy for assembling intricate structures at the micro- and nano-scale is to exploit self-assembly. By 
introducing colloidal building blocks into soft media such as fluid interfaces, nematic liquid crystals, or more complex 
mesophases [1-6], the resulting distortions of the medium can be used to fabricate more elaborate colloidal objects. 
Over the last ten years much effort has been devoted to achieving complex self-assembly by tailoring the shape of the 
elementary colloidal building blocks [5, 7]. Alternatively, the existence of mechanical instabilities in elastic materials 
suggests a simple but elegant method for guiding the creation of complex objects: a bottom-up assembly driven 
by elastic instabilities in simple, flexible structiires. Such a self-assembly might begin with simple building block 
configurations, such as straight filaments or flat sheets, vi^hich when coupled to another medium would fold or wrinkle 
to minimize the total energy of the system. If the system is well characterized theoretically, the dominant force 
balance can be tuned to yield desirable shapes. Recent examples include the use of adhesion and delamination [8-14] , 
and swelling and capillary interactions [15-17]. 

When particles are confined to a fluid interface, their interactions are mediated through the surface deformations 
[2, 18]. The deformations and thus the surface energy of the fluid may be decreased by a rearrangement of the 
particles on the surface, leading to effective capillary forces which can be attractive or repulsive depending on the 
geometry, density differences, and wetting properties of the particles. For two identical particles, the forces are always 
attractive. This so-called "Cheerios effect," and similar surface-mediated aggregation have been investigated in the 
context of vesiculation [19, 20], colloidal flocculation [2, 6], and millimetric ecology [21, 22]. While there have been 
many studies on filaments and polymers that interact via short-range forces and/or external flow flelds [23-31], the 
impact of long-range capillary forces on the shapes of flexible filaments lying on a fluid interface remains unexplored. 

In this paper we consider the problem of elastocapillary self-folding of a flexible micro-fllament. The model system, 
illustrated in Fig. 1, is that of a flexible filament confined to a fluid interface. The balance between long-range capillary 
forces, which are attractive between different parts of the filament, and bending stresses resisting deformation can 
lead to a novel buckling instability (as illustrated in Fig. lb), and at longer times to highly deformed folded patterns 
in the filament. Herein we show how to harness the shape instability to fabricate folded structures. The paper is 
organized as follows. In §11 we introduce the setup and discuss the energetics of the system. In §111 we pose the 
non-dimensional problem to solve for the shape dynamics in a viscous fluid. The linear stability analysis for the shape 
is presented in §IV, where we determine the functional dependence of the most unstable wavelength on both filament 
and fluid properties, as well as the effect of external tensile forces. We show that the threshold for self-buckling 
depends on a single dimensionless parameter quantifying the ratio between self-attraction and bending akin to the 
recently introduced elastocapillary number [15]. The complex, folded filament configurations reached far from the 
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FIG. 1: The model system considered in the paper, (a): A flexible filament (of radius a) floating at the interface between two 
fluids. Every part of the fllament experiences capillary attraction to every other part, while bending resists deformation. The 
competition between capillary forces and bending results in buckling and folding, (b): The most unstable wavelength depends 
on a dimensionless strength of self-attraction, Q. The onset of this instability can be tuned with the introduction of an added 
tension from an external source, F, which also decreases the most unstable wavelength (the shapes shown are obtained using 
the numerical computations). 
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II. ENERGY AND SCALING 



Consider two identical particles of size a located at the interface between two fluids, for example air and water. 
Provided the distance between the particles, R, is much smaller than the capillary length, their effective capillary 
interactions are described by the interaction energy ii^int = 70^ ln(i?/£c), where 7 is the interaction strength per unit 
area and the capillary length is defined by £c = \/o'/pg, with a the surface tension of the interface, p the difference in 
density between the two fluids, and g the acceleration due to gravity (£c « 2 mm for an air- water interface) [32] . The 
strength of the interaction comes from the competition between surface forces and an external force monopole [2]. 
This force can originate from gravity, as in the case of the "Cheerios effect" [18], or from electro- or magneto-static 
forces [4]; note that 7 contains additional information such as the particle geometry and wetting properties. 

The model system considered in this paper is illustrated in Fig. 1. We consider an inextensible clastic filament of 
length L < ic and bending modulus B. The filament bending energy is 
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where x(s, t) denotes the two-dimensional position vector of the filament on the fluid surface at time t. The arc-length 
coordinate parameterizing the fllament is denoted by s, subscripts denote partial derivatives, and T(s) is a Lagrange 
multiplier enforcing an inextensibility condition. For notational simplicity in the following we will commonly assume 
an implicit dependence upon time. The contribution to the energy due to capillary interactions may be written as 
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where i?(s, s') = |x(,s) — x(s')|, 7 is now the energy per area of the interaction between two points on the filament, 
and where we have defined /' = J^^'' + J^^a ' '^ith a the size of the short-length cutoff (a is typically of the order of 
the filament radius). 

Scaling distances by the filament length, L, the tension T{s) by B/L^, and energies by B/L, the total energy 
combining bending and interaction takes the simple form: 

i; = ^'|i|x,,|2+T(s)[|x,|2-l] + ^ j'\nR{s,s')ds'^ ds. (3) 

E depends on a single dimensionless parameter, = jL^/B, which quantifies the balance between capillarity self- 
attraction and bending, and is akin to the so-called elastocapillary number [15]. 

Next, we focus on deriving the equation of motion for a filament relaxing to equilibrium, which we will then exploit 
to analyze shape instabilities. 



III. EQUATIONS OF MOTION 



To perform detailed calculations of the filament dynamics, we consider variations of the total energy, E, with respect 
to the centerline position vector, x(s), and tension, T{s). The forces so derived must balance with the forces generated 
by the motion through the fluid, which we model through a local drag coefficient, C^t, an appropriate approach for 
dynamics in a very viscous fluid. The ambient fluid is assumed to be quiescent and of infinite extent. Scaling time 
upon the characteristic time scale r = L^/B, we arrive at a dimensionless equation for the motion of the filament, 

/•' xfs') - x(s) 

xt = -x,,,, +nj -S^J-^-^ds' + 5«(T(s)x,). (4) 

Anticipating our numerical approach, the local part of the integration in Eq. 4 may be handled analytically by adding 
and subtracting the singular part of the integrand. In so doing, we find 
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where c(s) = log[(l + S — s)/{s + S)], with S ~ a/L the dimensionless short wavelength cut-off. Rewriting the 
dimensionless filament position equation, we therefore have 

xt = -yissss + flc{s)xs{s) + nK[x{s')]{s) + ds (T(s)x,) . (7) 

Classically, the equation for the tension T(s) is derived using the relation Xg • x^^ = 0, and also the identities 

Xg * Xgg — 0, Xg • Xggg — Xgg • Xgg, (8) 
Xg * ^ssss 3Xgg • Xggg, Xg • Xggggg 4Xgg • Xgggg 3Xggg • Xggg, (9) 

yielding the differential equation 

(5gg - Ixggp) T{s) = - [4xgg • Xgggg + 3|xggg|2] - nc,{s) - nxg • J^K[x(s')](s). (10) 

Boundary conditions are derived by including terms to the energy functional to account for the work done in forcing 
the endpoints to move, namely, Ep = — Fq • x(0) — Fi • x(l), where Fq and Fi are (dimensionless) externally-applied 
forces at the s = and s = 1 endpoints of the filament, respectively. Assuming that the filament ends are free, hinged 
about a specified (possibly moving) point, or pulled by the forces defined above, the solvability conditions thus require 

Xgg(O) = xgg(l) = 0, (11) 
-Fo + Xggg(0)-T(0)x,,(0) = 0, (12) 
-Fi-Xggg(l)+T(l)xg(l) =0. (13) 
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FIG. 2: Instability of the filament shape: Exponential growth rate, a, computed numerically (symbols) and predicted analyti- 
cally (Eq. 25, lines), (a): With no external force, F = 0, we consider two interaction strengths Q. — 140 (just below the critical 
value fic = 144 for the instability), and Q = 2000. (b): With an external force F = 600, we consider the values Q = 2500 and 
Q, = 3200. 

The boundary conditions are then found to be 

x,,3(0) = (I - x,(0)x,(0)^) Fo, (14) 
x.,.(l) = - (I - x,(l)x,(l)^) Fi, (15) 
T(0)--Fo-x,(0), T(l)=Fi-x,(l). (16) 

To summarize, the partial differential equation describing the motion of a floating filament with either free or forced 
ends is given by Eq. 7 with the boundary conditions shown in Eqs. 11, 14, and 15. It is coupled to the equation for 
the tension, Eq. 10, with boundary conditions shown in Eq. 16. The numerical method developed to solve this system 
of equations is described in Appendix A. In the next section we exploit our equations of motion to describe the linear 
stability of the filament. 

IV. LINEAR STABILITY ANALYSIS 

A. Numerical investigation of growth rates 

Before we pursue an analytical description of filament instability, let us begin with a brief numerical investigation 
using the framework described in Appendix A. In order to consider the growth rate of a given wavenumber k, we seed a 
nearly straight filament with the perturbed dimensionless shape x(s) = {x{s),y{s)) = (s, oq sin(A;s)), with oq = 10"'*. 
Ten timesteps are then taken with a step-size At — 10~^, and the growth rate is computed as cr = log(aio/ao)/(10At), 
where aio = 2 sm{ks)y{s, lOAt) ds. Figure 2a shows (as symbols) the growth rates computed for two different 
interaction strengths when the filament ends are left free. For Q = 140, we find that there are no wavenumbers for 
which the growth rate is positive, indicating that the filament is linearly stable to transverse perturbations. Despite 
the attraction of every part of the filament to every other part, the filament's stiffness prevents the buckling and 
subsequent collapse when the body is perturbed from a straight conformation. Conversely, increasing the interaction 
strength to = 2000, the capillary attraction overcomes the bending elasticity and we observe that the filament is 
linearly unstable. A wide range of wavenumbers corresponds to positive growth rates. It is worth noting that the 
the critical value of the interaction strength above which capillary interactions destabilize the filament is much larger 
than 1, suggesting that it could not have been inferred solely from dimensional analysis. We also show in Fig. 2b the 
growth rates computed for two different interaction strengths, where an external pulling force has been introduced at 
the filament ends. For a given interaction strength the external force reduces the growth rate of every mode, and so 
decreases the value of the largest unstable mode. Meanwhile, setting F > introduces a negative slope in a for small 
k, and significantly changes the value of the smallest unstable mode. 
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B. Analytical predictions 

Having gained intuition from the numerical simulations described above, we now turn to an analytical investigation 
of filament instability. Our goal is to establish the dispersion relation governing the bending excitations of a floating 
elastic filament. 

One exact solution to the complete floating fllament system is the possibly unstable trivial solution, that of a 
straight filament lying along a director e, namely 

x*(s,t) = se-e- [Fi+Fo]i, (17) 

and 



T*{s) = e • [(sFi - (1 - s)Fo] + <| (1 - s) c(s) + log ( — ) + ^ log 



{l-s + S){s + S) 



S{1 + S) 



(18) 



Prom here, we are well-situated to perform a classical stability analysis based on perturbations of both x(s) and T{s). 
We assume that the endpoints are pulled upon by an external force of equal magnitude F, Fq = —Fe, Fi = Fe, and 

we set e = X for convenience of notation. 

To begin, we write x(s,t) = sx + ey{s,t)y + 0{e^), and let T = T* +T, with T* the equilibrium tension from 
Eq. 18. By tracking terms of size e and smaller, it may be shown that at leading order, the non-local contribution to 
the capillary force only gives rise to a transverse force contribution as 



K[x(s')](s) = 0(e')x + e 

Jo 



y{s') - y{s) jjs 



ds'y. (19) 



(s' — s)2 s' — s_ 

Upon the insertion of ^K[x(s')](s) = 0(e2)x + 0(e)y into Eq. 10, we deduce that T = T* + 0{e'^), such that we 
may dispose of tension perturbations T for the linear stability analysis. Similarly, we find that x • Xj = O(e^) from 
Eq. 7, so that the linear stability consideration may focus solely on transverse filament motions through the evolution 
oiy{s,t). 

In order to gain physical insight to the the nonlocal operator K given by Eq. 19, we assume that s is sufficiently 
far removed from the filament endpoints and look at its effect on plane waves of wavenumber k. By doing so, we find 
that the transverse component of K[x(s')](s) may be written as 

K[x(s')](s)-y«5(fc)y(s), (20) 

where 

g{k) = 2 [1 - cos(fc) + 7 + log(fc) - k Si{k) - Q(fc)] , (21) 
and with the sine and cosine integrals defined by 

Sm=[''^dt, a(^) = -r^rfi. (22) 
Jo t Jk 

The mathematical details of this derivation are provided in appendix B. Combining the results above, the transverse 
component of the filament position equation, Eq. 7, gives a linear relation for y(s), 

Vt = -yssss + ^9{k) y + [Oc(s) + T;] Vs + T*yss + 0{e^). (23) 

Upon inspection of T*, we are reminded that the local self-attraction term c(s) could just as well have been absorbed 
into the tension, and simplifying the above results in 

Vt = -Vsess + T*yse + O g{k) y. (24) 

We are thus left with the relatively simple over-damped Euler-Bernoulli beam equation, with a non-uniform tension, 
and with an additional term of magnitude fl originating from the non-local capillary attraction. In the following we 
discuss in detail the impact of this novel effective elasticity term on the dynamic response of the floating filament. 

To explore the response of the filament shape to small perturbations, we explicitly insert the ansatz y = 6**'''+'^*. 
While the finite-size effects must come into play for small wavenumbers and for longer times (trigonometric functions 
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do not satisfy the boundary conditions on Eq. 24, for example), this ansatz is expected to be accurate far from the 
end points of the filament. Upon insertion of the equilibrium tension T*(s), an inner product of Eq. 24 against e"**^* 
reveals the dispersion relation for the growth rate, 

a = -k^+^k^{n-2F) + ng{k). (25) 

How well does the analytical expression above capture the results of the full numerical simulation? Returning 
to Figs. 2a and 2b, the analytically predicted growth rates from Eq. 25 are shown as lines for each of the four 
cases considered. While discrepancies are expected due to end effects, which are not captured in our asymptotic 
consideration, in general we find good agreement between the numerical and analytical results. Therefore, we can use 
Eq. 25 with confidence to demonstrates the fundamental physics behind the elastocapillary instability, which we do 
below. 



C. Discussion of the elastocapillary instability 

The quartic term and quadratic terms in Eq. 25 are associated with the canonical Euler buckling. However, as 
expected, the capillary attraction leads to an additional compressive load (the flk^ term in Eq. 25). Any attractive 
short-range coupling along the filament is expected to renormalize the tension term in a similar fashion. We now 
discuss the consequence of the more complex contribution: flg{k), to the dispersion relation. We first stress that 
^(0) = 0, and g{k) < for all k > 0, which means that the non-local part of the capillary coupling acts to stabilize 
straight conformations. 

In addition, given that O > and F > 0, the predicted growth rate a goes to zero as fc ^ 0, to — oo as A; ^ oo, 
and may or may not take on positive values for finite fc, depending on the values of fi and F. For example, Fig. 3a 
shows the growth rate from Eq. 25 as a function of wavenumbcr k for three values of the interaction strength, il, with 
an external pulling force F = 200. For fl less than a critical value, flc = 1138 in this case, all modes are predicted 
to decay in time, and the filament is predicted to be stable. For Q > fig, we observe the appearance of a band of 
unstable modes and the filament becomes unstable. The critical interaction strength for instability, ilc, is shown as a 
function of the value of F in Fig. 3b. It takes the value = 144 when F = and increases with increasing F. Both 
the values predicted by the analytical theory and those found by the full numerical simulation are shown; here again 
we see that the analysis is sufficient to capture the linear instability of the filament dynamics. 

To develop physical insight into this behavior, we look at at the asymptotic form of (7{k) in the limit of long 
wavelengths. We find that a ~ —Fk^ + (f2/144 — l)fc** for fc <C 1. In this regime, the capillary-induced compressive 
load vanishes; the self-attraction does not renormalize the tension of the filament in the long wavelength limit. 
Nevertheless, the capillary coupling softens the filament by effectively reducing its bending stiffness: the pre-factor of 
the quartic term decreases linearly with and changes its sign at flc = 144. Above this critical value the bending 
excitations along the filament are amplified. As noted in the previous section, this unexpectedly high value of ilc 
could not have been anticipated on the basis of dimensional arguments alone. Furthermore, as da/dk < as A; — >■ 
for F > the smallest unstable mode, if there is one, must take on a finite value when an external pulling force is 
applied. 

Beyond this prediction of the destabilization threshold, we can also accurately predict the value of the most unstable 
mode, which chiefly sets the ultimate folded shape of the filament as we will show more precisely in the next section. 
For a given force F and interaction strength f7 (such that > fJc), the most unstable mode k* is recovered upon 
setting da/dk to zero and solving for k. The resulting value is presented in Fig. 3c as a function of both Q. and F. 
Predictably, the most unstable mode increases with increasing self-attraction (increasing il), and decreases as the 
filament ends are pulled upon with greater intensity (increasing F). The introduction of an external force renders 
the filament stable for small interaction strengths f2, and so intuitively prolongs the onset of buckling. Increasing 
Q. and/or decreasing F is therefore expected to yield a more finely folded filament shape once the instability fully 
develops. 

In the specific case of a free filament, F = 0, Fig. 3d shows the value of most unstable wavenumber as a function 
of ri. The shaded region indicates the band of unstable modes for a given interaction strength Vl. Again, these 
variations can be qualitatively understood by inspecting the asymptotic behavior of a{k). Indeed, for fc 1 we find 
(j{k) ~ — fc* + f2 [fc^ — Trfc -F 2(log(fc) + 7 + 1)] , with 7 the Euler-Mascheroni constant. This latter approximation of 
(j(fc) is accurate to within 2% of the exact value for wavenumbers as small as A = tt, and it can therefore replace 
(j{k) generally for all but the very longest perturbation wavelengths. For k 1, the problem boils down to the 
canonical Euler Buckling induced by the longitudinal self-attraction. Hence, balancing the cubic and linear terms in 
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FIG. 3: (a): Growth rate of the instabihty, a, as a function of wavenumber, k, for three values of the interaction strength, Jl, 
with external force F — 200 (the value of the critical interaction strength is Jlc = 1138). (b): Critical interaction strength, fie, 
below which there are no unstable modes as a function of the external force, F, from both numerical simulations (symbols) 
and analytical prediction (line), (c): Map of the most unstable mode, k* /tt, analytically predicted (from Eq. 25) as a function 
of and F. (d) Setting F = Q (free filament ends), the most unstable wavenumber, k* , is surrounded by an increasingly wide 
band of unstable modes as Q. increases (shaded region); analytical estimates for the most unstable and the highest unstable 
mode are shown as dashed lines. 



the equation da/dk — gives the usual approximation to the most unstable mode [33]: 

k* « ^Vf^ - 2F. (26) 

The largest unstable mode may also be predicted, again assuming /c 1, by setting a ~ — fc^ + (ri/2 — F)k'^ = 0, 
yielding k^ax ~ — 2F / y/2 = y/2 k* . These approximations are shown as dashed lines in Fig. 3d. Both overestimate 
their true values in the present study but capture the numerical resuhs qualitatively. The reason for this discrepancy 
is due to the range of relevant wavenumbers involved in the elastocapillary instability. Destabihzation of the bending 
modes occurs at small or intermediate values of the wavenumber fc; the non-local corrections to the Euler-buckling 
problem {flg{k) terms in Eq. 25) thus cannot be overlooked and strongly contribute to the dynamic response of 
floating microfilaments. 

In summary, we have shown in this section through a linear stability analysis that a self-attracting filament displays 
a novel mechanical instability. We have explored the structure of this instability, and its dependence on both the 
strength of the self-attraction and the stabilizing external force. 



V. BEYOND LINEAR-STABILITY: SELF-FOLDING OF A FLOATING FILAMENT 



The linear stability analysis of the previous section yields some insight into how folding will begin. In this section we 
investigate long-term, large-amplitude folding behavior. The path that a filament takes towards any final configuration 
depends on the initial conditions, the material parameters of the system (specified completely by fi), and the external 
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FIG. 4: Long-time dynamics of filament shape, (a): With Q = 2500 and F = 0, two initial perturbations with different initial 
wavenumbers (fc = vr, top and k = 7tt, bottom) fold into similar structures but at different times, (b): Adding an external 
force, F, modifies the folding shape of the filament; with fl = 2500 and F — 500, the final folded shape has fewer loops than 
in the F = case due to the straightening influence of the added tension. Meanwhile, for F = 800 a secondary instability 
develops along the axial direction, leading to a final left-right asymmetric shape, (c): With SI = 3500, a highly deformed initial 
hairpin configuration undergoes a cascade of self-folding events, (d): Folding and unfolding of a floating filament is shown for 
n = 3500; the fllament self-folds with F = until t = 9 x 10~^ at which time a pulling force F = 3500 is introduced and the 
filament unzips, albeit through a different path through shape space. 

tension on the filament (F). In this section we compute the dynamics from the equations of motion described in §111 
using the numerical method described in Appendix B. 

For large values of the interaction parameter the capillary attraction overcomes the relatively weak bending re- 
sistance, and the filament displays a tendency to "accordion": it collapses with a short wavelength relative to the 
filament length. The dynamics of such filament self-folding are shown in Fig. 4a, where we plot the value of the 
total energy E (Eq. 3) as a function of time, fixing — 2500 and F — for two different initial filament shapes. 
The topmost curve (blue online) shows the path taken of a filament that is seeded with an initial transverse position 
y{s, 0) — 0.05 sin(7rs) as it travels down to a state of lower energy. At short times the shape is predominately composed 
of a number of anti-nodes corresponding to the fastest growing wave number, predicted by maximizing Eq. 25 to be 
k* = 7.47r. These anti-nodes grow to form loops, and eventually the filament approaches a regime of near-contact. 
In this vicinity there is a slight discontinuity in the slope of the energy, as a result of introducing near-field repulsive 
physics required to regularize the diverging capillary interactions (see Appendix A). The lower curve in Fig. 4a (red 
online) shows the path taken by an identical filament seeded with a different initial condition, y{s,0) = 0.05 sin(77rs). 
Here the filament eventually achieves the same conformation as in the previously considered case, dictated by the 
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fastest growing linear mode. However, due to its initial shape being closer to the fastest growing linearly unstable 
mode at the outset, this filament begins its descent into lower energy conformations earlier. Interestingly, this result 
suggests that the elastocapillary buckling provides a robust means to achieve complex foldings of microfilaments. 
Indeed, the gross features of the ultimate shapes can be anticipated from a simple linear stability analysis. 

By introducing an external tension, _F > 0, the possibilities for folding become more versatile but more complex. 
In Fig. 4b we show three different filaments, all with initial positions given by y{s,0) = 0.05sin(7rs + 0.01s) (a shape 
with a slight left-right asymmetry), but pulled upon with three different magnitudes of external forcing (F = 0, 500, 
and 800). With the addition of this external tension the most unstable mode structure is changed, and a different 
final state is selected. With no external tension {F — 0), the filament accordions tightly as in Fig. 4a. In the second 
case, with F = 500, the fastest growing unstable wavenumber is shifted to a smaller value, as predicted by the linear 
theory, leading to a smaller number of loops in an equilibrium conformation. Finally, with F = 800 we observe a 
secondary instability growing along the axial direction, perpendicular to the transverse instability at the center of the 
previous section. 

In fact, when a filament is sufficiently long relative to the structures generated by the initial instability, it can 
undergo a cascade of instabilities. This is illustrated in Fig. 4c, where we seed the filament initially with the shape 
of a hairpin. In the absence of external forcing, the hairpin shape is seen to first "zip up" into a racquet shape, 
reminiscent of the ones observed in wetted fibers [15], and in semi-flexible polymer systems [23, 24]. Over time it 
displays a cascade of buckling and folding. By pulling on the filament ends with an external force, F = 500, similar 
multiply folded states are obtained, yet the final shape is markedly different from the one observed in the force-free 
case. 

Finally, in Fig. 4d we illustrate the unfolding of a self-folded structure. The filament is set initially with the shape 
y{s,0) = 0.05sin(77rs), and we fix = 3500 and F = 0. The shape evolves in a similar fashion as in the previously 
described examples (blue solid line). At t = 9 x 10~^, we introduce external forcing on the filament ends with 
magnitude F = 3500 to unfold (black dashed line). We observe an unzipping of the folded structure, with jumps in 
the filament energy indicative of each loop being pulled apart in succession. Rather than a reversible shape change 
we clearly observe unfolding along a different path. 

VI. CONCLUSION 

In this paper we have described the elastocapillary instability of flexible filaments confined to a fluid surface. 
Buckling and folding result from the interplay between self-attraction through long-range capillary forces and elastic 
bending rigidity. As a result of the competition between these mechanical forces, the instability only appears beyond a 
large, critical ratio of the capillary self-attraction to the filament's bending modulus. The most unstable wavenumbers 
are typically located at an intermediate range, thereby yielding the spontaneous formation of accordions folds. Due 
to the the long-range nature of the capillary interactions, both the wavelength and the growth rate of this instability 
differ from the conventional Eulcr-buckling results, which would be observed for short-range self-attraction. 

After having analyzed the onset of filament shape instability through a linear theory, numerical simulation was used 
to show a number of possible conformational transitions for the self- folding filament. In particular, introducing an 
external pulling force at the filament ends changes the unstable linear modes, and modifies significantly the long-time 
behavior. For a filament sufficiently long relative to the shapes driven by self-attraction there can be a cascade of 
folding instabilities, each acting perpendicularly to the previous fold. While the linear analysis is only valid for short 
times and nearly straight filaments, our simulations demonstrate a wide variety of accessible conformations, from 
rackets and sewing needles to accordion shapes. 

The numerical method described here could be augmented to study many related problems. Future work might 
include computing the relative energy scales for each equilibrium shape; with sufficient statistics, the ensemble of 
different shapes that could be expected should our system be coupled to a thermal bath might be elucidated. The 
deterministic simulations could also be used to determine the statistical significance of certain shapes in crumpled 
thin sheets or fibers [34, 35]. 

Our original experimental motivation was provided by buoyant particle interactions, where the dominant competi- 
tion is between elasticity and the density of the particles. In which practical experimental setup could the buckling 
instability be observed? We recall that the critical elastocapillary number is given by f2 = ■yL^/B. The instability 
is cixpected to be observable for filaments both soft (small value of B) and dense (large value of 7). Consider a 
filament composed of microns-sized spherical colloids linked by double-stranded DNAs as in Ref. [36] . From capillary 
interactions we expect 7 ~ F'^/a where F is the force per unit length acting on the filament and a the surface tension 
of the fluid interface. From gravity we get F ^ Apa^g and thus 7 ^ Ap^a^g'^ /<j. For colloids of radius a = 10 /im 
and Ap = lO'^kg/m^ and an interface with a ~ 10~^ N/m we would get 7 ~ 10~^° N/m. Using filaments of length 
L 1 mm, right below the capillary length, and with a DNA linker flexibility of B 10~^^ J/m [36], we would obtain 
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FIG. 5: A near-contact situation, s* denotes the point on the filament closest to the station s. 



an elastocapillary number or order f2 ~ 10'^, and thus the instabiUty could be observable. Obviously, beyond the case 
of gravitational forcing, any external force acting on the filament can serve to deform the interface and drive similar 
elastocapillary self-folding. Electric or magnetic fields acting on self-assembled filaments of metallic nano-particlcs or 
super-paramagnetic nano filaments could be two experimental examples. 
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Appendix A: Numerical solution 

The numerical approach taken in this manuscript is a variation on that described by Tornberg & Shelley [37]. 
The complications present here are the stiffness introduced by the elastic restorative term cx x^sss, and the nonlocal 
integral operator K[x(s')]. The approach will be to treat the stiffest terms implicitly in time, while treating the 
remaining terms explicitly for a fast but stable scheme. 

Time is discretized uniformly, with the n*'* timestep denoted by tn, and writing (for instance) the tension at the 
n*'* timestep in shorthand as T{s,tn) = T"(s). At each time i„ the tension is first solved via 

(a,, - T«(5) = -A [4x';, . x^,,, -I- 3\^:,f] - c,(s)- (Al) 

x^. J^K[x"(s')](s)+/3(l-x,.x,). 

We have included a correction term to remove numerical errors in the determination of the arc-length, again following 
[37]. To solve this equation we define a uniform grid s = i/M, for i = 0, 1, M. Spatial derivatives computed using 
second-order accurate finite difference formulae. Since the position x" at time tn is known, the equation above is 
inverted to recover T„(s) without difficulty. 

Solving the filament position equation is the next step, where we update via 

^ (3x"+i - 4x" + x"-i) = -x«+i + nc{s){2^- - xri)+ (A2) 

17K[2x"(s') - x"-i(s')](s) + ds (2T"x^ - r'-^x^^^) . 

For the first time step the above is replaced by a first order Euler step. In this case, we need only build the following 
matrix operator and to invert it once: M = 3/(2At)I -I- (along with boundary conditions). The inverted operator 
may be applied to known quantities at every subsequent time step, providing for a very efficient numerical scheme. 
We have subtracted off the singular part inside the integrand of the interaction forcing so that the remainder of the 
expression is finite, as described in the main text. 

A challenging numerical problem is introduced when the filament comes into near contact with itself. Or worse, 
when the filament exhibits a large number of loops so that there are many near-contact interactions. The situation is 
illustrated in Fig. 5, where two different parts of the filament are shown nearing each other; the centerline is denoted 
by dotted lines, and the blue solid lines indicate the thickness of the filament. 

Our intention is to prevent filament crossing, but to allow for tangential motions in the near-contact state. For 
each station in arc-length on the filament s, we first determine the closest point on the filament s* that is at least a 
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distance 4e from s (so as to neglect neighboring points on the filament), 

s*(s) =min{|x(s') -x(s)| ,|s'-s|>4e}. (A3) 

s' 

A local coordinate system is then defined by 0{s,s*) = (x(s*) — x(s))/|x(s*) — x(s)| and fj = v-^ . The forces in 
the filament equation are then decomposed locally into their u and ry components. The tangential component of the 
filament force balance is left unchanged, but the normal component is adjusted as follows 

i>{s, s*) ■ (xt + :x.ssss) = (A4) 

min{i>(s,.*).G,P(.,0-G[l-( ^^^^^^;^^,J ^]} 

where G denotes the right hand side of Eq. 7 (minus the stifFcst term, x,,,,), and p > (we use p = 12 for all the 
results presented here) . The multiplication factor is such that there is no normal component of force (outside of the 
elastic part) when the filament centerline is 2e away from the closest point on the other part of the filament centerline. 
The minimum of the two quantities is chosen so that if the fiuid and tension forces wish to separate the filaments 
from each other they can, but if the fluid and tension forces are pushing the filaments to an overlapped state the sign 
of this force is reversed. This treatment of the near-contact physics appears to allow the filament to fold up tightly 
without numerical difficulties, without significantly impacting the underlying physics of capillary attraction. When 
the filament comes into near-contact, the timestep size is decreased and the spatial resolution is increased adaptively. 

For the linear stability analysis we set M = 2000, At = 10~^, /3 = 5 * 10^. For the long-time dynamics and 
equilibrium shapes it was sufficient to set M = 600, A.t = 10~*, /3 = 5 * 10^, and p = 12, with adaptive refinement of 
the grid and decrease in timestep size as necessary when the filament approached itself as described. 



Appendix B: Nonlocal integration 



We begin by showing, for fc ^ 1 and s sufficiently far removed from the filament endpoints, that 



(Bl) 



for y{s) oc sin(fcs). For fc ^ 1, we are led to write the integration over a new inner variable s' = s + £,/k, so that 

fk(l-s) 



ms)] = k [ { ^^^(fe^O - ^ } d^, (B2) 

J — ks S s 



and for s well distanced from both endpoints, we may approximate the above as 



=.sin(M r S^d^+.cos(M r ^'-^^-^ 

J — OO > J — oo 



= —Trksm{ks) = —nky{s). (B3) 

Though we do not reproduce it here, a complicated analytical expression for K[y{s)] in Eq. Bl may in fact be derived 
for arbitrary wavenumbcr k. While it is not generally the case that i^[y(s)] is diagonalized using a trigonometric basis 
(as it is in the k^ 1 approximation above), we still note that 



/" 



r-l 

e-iksK[e''''']ds = g{k), (B4) 

/o 

gik) = 2[1 - cos(fc) +-f + log(fc) - k S,ik) - Ciik)], (B5) 



with the sin and cosine integrals defined in Eq. 22. For fc ^ 1 we have g{k) ~ —k^/2, while for fc 1 we have 
g{k) ~ —nk + 2(log(fc) +7 + 1), with 7 the Euler-Mascheroni constant. This latter approximation of g{k) is accurate 
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to within 2% for wavenumbers as small as k = n. 
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